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1. Introduction 

■ The physics of matter at subnuclear densities p G [10^^ — 10^''] g cm ^ is of great 
interest for the astrophysics of compact objects. The "hot" stage of evolution of 
matter, in which temperatures are in the range of tens of MeV, is associated with the 
dynamics of supernova explosions. Knowledge of the equation of state, composition, 

I and weak-interaction processes are of prime importance for an understanding the 

mechanism of explosion, the formation of neutrino spectra at the neutrinosphere, 
and the elemental abundances of the low-density matter in the supernova winds 
that are prerequisite for the onset of r-process nucleosynthesis. Days to weeks after 
the supernova explosion subnuclear matter has become "cold" , with temperatures 
T < 0.1 MeV. Moreover, the properties of the subnuclear matter forming the crust of 
a neutron star are of fundamental importance for the entire spectrum of observable 
manifestations of pulsars, ranging for example from superfluid rotation dynamics 
to magnetic field evolution to neutrino cooling. 

Nuclear matter at subnuclear densities is a strongly correlated system in which 
the relevant degrees of freedom are well established and the interactions are con- 
strained by experiment. The challenge lies in the many-body treatment of this sys- 
tem where macroscopic quantum phenomena such as Bose-Einstein condensation of 
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deuterons and alpha particles exist as well as the BCS pairing in neutron matter. 
Our aim here is to describe some of the many-body methods for dealing with such 
correlated states of matter. We will pay less attention to the physical setting and 
implications of the results; the reader concerned with these issues is referred to the 
original literature cited among the references. 

2. Bose-Einstein condensation: a lattice Monte-Carlo perspective 

In this section we describe an approach to interacting Bose systems which is valid 
in the vicinity of the critical temperature of Bose-Einstein condensation (BEG). 
The method was put forward in the context of dilute gases interacting via repulsive 
two-body forces^ and has since been reformulated for a strongly correlated system 
interacting with attractive two-body and repulsive three-body forces.'^ The method 
has been applied to Bose condensation of alpha particles in infinite matter. (Alter- 
native studies are based on hypernetted-chain summations."*) Consider a uniform, 
non-relativistic system of identical bosons described by the Hamiltonian 

^VV;t(x)VV'(x) - m|V'(x)P + 52|V^(x)|^ -t- 93m^)\' , 
Zm 

where m is the alpha-particle mass, /i is the chemical potential, and ijj is the boson 
field. Below, we shall implement lattice regularization. The theory defined by Eq. ([T]) 
can be mapped onto an effective scalar field theory within the finite-temperature 
Matsubara formalism. Consider the fields tp and tp^ as periodic functions of the 
imaginary time t S [—/?,/?], where (3 = 1/T is the inverse temperature. Next, 
decompose the fields into discrete Fourier series 

oo oo 

V(x,c.,)= J2 e'"^^V(x,T) = V^o(x)+ J2 e'""'V'(x,T), (1) 

1^= — oo i/=~oo, i/T^O 

where the Fourier frequencies are the bosonic Matsubara modes lUi, ~ ^ixivT 
(with V taking integer values). The Matsubara Green's function is given by 
G'^^(a;i,,x) = [icjy — (2r7i)^*V^ + Here the chemical potential may include any 
contribution from the momentum- and energy-independent part of the self-energy; 
we also assume that any momentum and energy dependent parts are absorbed in 
the mass and the wave-function renormalizations, respectively. Since /i — > near Tc, 
the characteristic scales of spatial variations of the Green's function with non-zero 
Matsubara frequencies are I — {2mijjy)~'^/'^ , which are of the order of the thermal 
wave-length A = {2n/mTy^^. The contribution of the non-zero Matsubara modes 
to the sum in Eq. ^ will be neglected since we are interested in scales L I, which 
are characterized only by the zero- frequency modes. In terms of new real scalar fields 
(f>i and 4>2 defined via the relations "00 = v{4>i + *</'2) and ipo — vi'Pi ~ where 
rj = m/h?(3, the continuum action of the theory is given by 

S{cj^)= j d3^|i^[a.0(x)]V^0(x)2- J[0(x)2]Vg[</.(x)2]n, (2) 
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where (fP = cfif + (f>2, r = —2f3iJ,ri^, u = Al/3g2ri* and w = &.(3g3rf. The action 
([2]) describes a classical 0(2) symmetric scalar 0^ field theory in three spatial di- 
mensions (3D). The positive sextic interaction guarantees that the energy is bound 
from below, which would not otherwise be the case because of the negative sign of 
the quartic term describing the attractive two-body interactions. The characteristic 
length scale of the theory is set by the parameter u, which has the dimension of 
inverse length; the dimensionless parameter of the lattice theory is ua^, where is 
the lattice spacing. The thermodynamic functions of the model are obtained from 
the partition function 

Z = j [d</)(x)]exp[-5(^)]. (3) 

For example, the expectation value of the particle number density is given by Uq = 
= {/3V)~^dlnZ/dfi, where V is the volume. The continuum theory is now 
discretized on a lattice by replacing the integrations over spatial coordinates by a 
summation over lattice sites. The discretized version of the continuum action ^ is 

Sl (-^^ = E I E W<^i + <^^)M^f + A [1 + M^f] ' A + C [M^)^] ' I , 

the v summation being carried out over unit vectors in three spatial directions 
(nearest neighbor summation). The hopping parameter n and the two- and three- 
body coupling constants A and C arc related to the parameters of the continuum 
action through a|r = (1 — 2X)/k — 6, X — and ( — wk^/90. The lattice 

and continuum fields are related by (/'l(x) = (2K/ai,)^/^0(x). The components of 
the spatial vector discretized at integral multiples of the lattice spacing gl- 

Xfi = 0, ol, . . . (Li/ — l)ai,. For a simple cubic lattice in 3D with periodic boundary 
conditions imposed on the field variable, one has 4)l{x + a^L) = 4>l{x) (N.B. for a 
box of length L there are (real) variables within the volume {Lgl)^). 

In Ref. 3, the field configurations on the lattice were evolved using a combination 
of the heatbath and local Metropolis algorithms, executing 10^ — 10^ equilibration 
sweeps for lattices sizes from 8'^ to 64'^. Once the field values on the lattice were 
determined, these were transformed into their counterparts in the continuum the- 
ory to obtain the statistical average value (H) of the Hamiltonian, i.e., the grand 
canonical (thermodynamical) potential 17 as a function of density. The critical tem- 
perature Tc for Bose-Einstein condensation can be obtained from the simulations.'^ 
In practice one obtains the density n{T), rather than T(n) directly, at constant fu- 
gacity z 1. Working at small fugacity log z = —0.1, the density of the system at 
constant small chemical potential can be computed and the associated temperature 
is then identified with Tc- 

3. From pair condensation to three-body bound states 

The few-body bound states are of interest in the "hot" stage of compact stars to the 
extent that they can provide a very efhcicnt source of opacity for neutrinos propa- 
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gating through matter. This situation can be compared to that of radiative transfer 
in ordinary stars, where the photoabsorption on the weakly-bound negative ion of 
hydrogen (H^) largely determines the opacity. If the system is isospin-symmetric, 
the pairing occurs in the '^Si Di partial- wave channel. There is a two-body 
bound state in this channel in free space - the deuteron; hence the BCS to BEC 
crossover arises in this new context. Following up on the conjecture of Nozieres 
and Schmitt-Rink,^ one may attempt to describe this crossover within mean-field 
BCS theory. The central numerical problem then reduces to solution of the gap 
equation 

A;(p) = _ / y v-iSDi^ ,^ ^1'{P') [1 -2f{E{p'))] , (4) 

where D^{k) = (3/87r)[AQ(fc) -I- A|(fc)] is the angle-averaged neutron-proton gap 
function, V^^^'^{p,p') is the interaction in the "^Si —^Di channel, E{p) is the quasi- 
particle spectrum, and / is the Fermi distribution function. The chemical potential 
is determined self-consistently from the gap equation (j4|) and the expression for the 
density. 

We now outline an algorithm for numerical solution of the gap equation, which 
can be applied to arbitrary potentials that are attractive at large separations.^ The 
method is effective in dealing with the hard core (short-range repulsion) in nuclear 
potentials and could be useful for other systems featuring short-range repulsive 
interactions. The starting point is the gap equation with an ultraviolet momentum 
cutoff A ^ Ap, where Ap is of the order of the natural (soft) cutoff of the potential. 
Successive iterations, which generate approximant Aj to the gap function from 
approximant A*^*^^^ (i = 1, 2, . . .), are determined from 

A<-)(P,A)= r^y3S.i(^^^,) A(-^)(y,A) [,_,^(^(^))]. 

The process is initialized by first solving Eq. ([4]) for D{pp), where pp is the 
Fermi momentum, assuming the gap function to be a constant. The initial ap- 
proximant for the momentum-dependent gap function is then taken as A''^'^^(p) = 
V{pftP)D{pp). Two iteration loops are implemented at given chemical potential. 
An internal loop operates at fixed A and solves the gap equation ^ iteratively 
for « = 1, 2, . . . . An external loop increments the cutoff A until dA{p, A)/dA 0. 
The finite range of the potential guarantees that the external loop converges once 
the entire momentum range spanned by the potential is covered. Thus, choosing 
the starting A small enough, we execute the internal loop by inserting A^*~^^(p) in 
the right-hand side of Eq. ([5]) to obtain a new Ai(p) on the left-hand side, which 
in turn is re-inserted in the right-hand side. This procedure converges rapidly to a 
momentum-dependent solution for the gap equation for A{p, Aj)6{Aj —p), where 9 
is the step function and the integer j counts the iterations in the external loop. 

For the next iteration, the cutoff is incremented to Aj = Aj_i + (5A, where 
5 A <C Aj, and the internal loop is iterated until convergence is reached. The two- 
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loop procedure is continued until Aj > Ap, after which the iteration is stopped, 
a final result for A{p) independent of the cutoff having been achieved. Once this 
process is complete, the chemical potential must be updated via the equation for 

the density. Accordingly, a third loop of iterations seeks convergence between the 
the output gap function and the chemical potential, such that the starting density 
is reproduced. 

The BCS-BEC crossover^ in nuclear systems has been studied as a function of 
density and asymmetry in the population of isospin states (proton- neutron asymme- 
try). A number of interesting features are revealed. In the extreme low-density 
limit, the chemical potential changes its sign and tends to —1.1 MeV, half the 
deutcron binding energy. Thus, twice the chemical potential plays the role of the 
eigenvalue in the Schrodinger equation for two-body bound states in this, the BEG 
limit. The pair function is very broad in this low-density regime, indicating that the 
deuterons are well localized in space; conversely it is peaked in the BCS limit at 
high density where the Cooper pairs are correlated over large distances. In the case 
of asymmetric systems, the density distribution of the minority particles (protons) 
has a zero-occupation (blocking) region that is localized around their Fermi-surface. 
Upon crossover to the BEC side, the blocking region becomes wider and moves to- 
ward lower momenta. Eventually there is a topological change in the Fermi surface: 
the states are occupied starting at some finite momentum and are empty below 
that point. The Nozieres-Schmitt-Rink conjecture^ of a smooth crossover from the 
BCS to the BEC limit does not hold in general for asymmetric systems; instead, 
phases with broken space symmetries intervene within a certain range of population 
asymmetries (see Ref. 9 and references therein). 

The three-body bound states can be computed from the three-body scattering 
matrix, which is written as T = T^^) -f-T^^^ +T'^^\ with components defined (using 
operator notation for compactness) as 

03 

' fl — ei — 62 — €3 + ir] 

These are nonsingular type II Fredholm integral equations; the operator Qs in mo- 
mentum representation is given by (53(ki, k2, ks) = [1 — /(ki)] [1 — /(k2)] [1 — /(ks)] — 
/(ki)/(k2)/(k3); and ei(ki) are the quasiparticle spectra. The momentum space 
for the three-body problem is conveniently spanned by the Jacobi four-momenta 
K = ki + kj + kk, Pij = {ki — kj)/2, and qk = {h + kj)/3 — 2kk/3. The Tij-matrices 
are essentially the two-body scattering amplitudes, embedded in the Hilbert space 
of three-body states. In the momentum representation they are determined from 

(7) 

where Q2{q,p) = (1 - /(q/2 + p) - /(q/2 - p)) and €±{q,p) = (e(q/2 ± p)) are 
averaged over the angle between the vectors q and p. Compared to the free-space 
problem, the three-body equations in the background medium now include two- 



= + %j- ^ — (rW + T(^)) . (6) 

■' il - 61 - eo - eri + 111 \ ) 
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Fig. 1. Dependence of the two-body (E^) and three-body (i?t) binding energies on inverse 
temperature, for fixed values of the ratio / = no/n, where n is the baryon density and no = 0.16 
fm~^ is saturation density of nuclear matter. For asymptotically large temperature, Eii{oo) = 
-2.23 MeV and Et{oo) = -7.53 MeV. The ratio Et{P) / Ed{P) is a universal constant independent 
of temperature.* 



and three-body propagators that account for (i) the suppression of the phase-space 
available for scattering in intermediate two-body states, encoded in the functions 
Q2, (ii) the phase-space occupation for the intermediate three-body states, encoded 
in the function Qs, and (iii) renormalization of the single-particle energies t{p). 
For small temperatures the quantum degeneracy is large and the first two factors 
significantly reduce the binding energy of a three-body bound state; at a critical 
temperature T^z corresponding to Et{(3) — 0, the bound state enters the continuum. 

This behavior is illustrated in Fig. [U which shows the temperature dependence 
of the two- and three-body bound-state energies in dilute nuclear matter for several 
values of the density of the environment. In analogy to the behavior of the in- 
medium three-body bound state, the binding energy of the two-body bound state 
enters the continumn at a critical temperature Tc2, corresponding to the condition 
Ed{/3) = 0. Our solutions exhibit a remarkable feature: the ratio rj = Et{(3) / Ed{/3) 
is a constant independent of temperature. For the chosen potentials, the asymptotic 
free-space values of the binding energies are Et{0) — —7.53 MeV and Ed{0) — —2.23 
MeV; hence 77 = 3.38. An alternative definition of the critical temperature for 
trimer extinction is Et{P'^^) = Ed{(3). This definition takes into account the break- 
up channel t — > d + n oi the three-body bound state into the two-body bound 
state d and a nucleon n. The difference between the two definitions is numerically 
insignificant. 

Fig. [2] depicts the normalized three-body bound-state wave function for three 
representative temperatures, as a function of the Jacobi momenta p and q. As 
the temperature drops, the wave function becomes increasingly localized around 
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Fig. 2. Wave function of the three-body bound state as a function of the Jacobi momenta p and 
q defined in the text, for / = no/n = 60 and temperatures T = 60 (left panel) and 6.6 MeV (right 
panel). 



the origin in momentum space. Correspondingly, the radius of the bound state 
increases in r-space, eventually tending to infinity at the transition. The wave- 
function oscillates near the transition temperature (right panel of Fig. [2]). This 
oscillatory behavior is a precursor of the transition to the continuum, which in the 
absence of a trimer-trimer interaction is characterized by plane-wave states. 



4. The weak interaction in cold subnuclear matter 

Non-nucleonic channels of cooling that operate in the crusts of neutron stars are 
electron neutrino bremsstrahlung off nuclei and plasmon decay:^" e + {A^ Z) -^■ 
e -|- {A^ Z) + ly + D and plasmon +1/ + D. Above the critical temperature Tc for 
neutron superfluidity, the neutrons that occupy continuum states (i.e., those not 
bound in clusters) emit neutrinos of all flavors / via the bremsstrahlung process^^ 
n + n^n + n + i/f + Df. At T<Tc the latter process is suppressed exponentially by 
exp(— 2A/T), where A is the gap in the quasiparticle spectrum. The superfluid na- 
ture of the matter allows for a neutrino-generating reaction (known as pair-breaking 
and recombination), whose rate scales like and thus is specific to the superfluid 
(i.e., vanishes as A — > 0). The rate of the process is given by the polarization tensor 
of superfluid matter. 

A systematic diagrammatic method to compute the reaction rates is based on 
the kinetic equation for neutrino transport, formulated in terms of real-time Green's 
functions. "'^^ The corresponding Boltzmann equation is 

[dt + d,u;M)d.]Mci,x) [ ^Tr[n<{q,x)S>{q,x)-n>{q,x)S<{q,x)], 

Jo '^'^ 

where q = (gojQ) is the four momentum, S'^'^(q,a;) are the neutrino propagators, 
and f2^'^(g,x) are their self-energies. In second Born approximation with respect 
to the weak vertices , the latter are given in terms of the polarization tensor (s) 
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Fig. 3. The sum of polarization tensors that contribute to the neutrino emission rate. The con- 
tributions form n('')(g) and Il^'^\q) vanish at the one-loop approximation. 

n^A^ {Qi ) x) of ambient matter as 

(8) 

The "greater" and "lesser" signs refer to the ordering of two-point functions along 
the Schwingcr contour in the standard way. 

It is the total loss of energy in neutrinos per unit time and unit volume, i.e., 
the emissivity, that is of interest for the astrophysics of compact stars. This quan- 
tity is obtained by integrating the first moment of Boltzmann equation. For the 
bremsstrahlung of neutrinos and anti-neutrinos of given flavor it is expressed as 

q2 I- 

f^vv = — ^ / ^^"9^(91 + 92 - 9)9o5f(9o)A''^(gi,92)Imn^c(9). (9) 

where G is the weak coupling constant, q is the four-momentum transfer, g{qo) = 
[ex.p{qo/T) — is the Bose distribution function, n^^(g) is the retarded polar- 
ization tensor, and A'^^{qi,q2) = Tr [7'^(1 — 76) ^i7''(l — 76) ^2]- Sums over the 
neutrino momenta (71.2 indicate integration over the invariant phase-space volume. 
The central problem of the theory is to compute the polarization tensor of the cold 
subnuclear matter. Initially, calculations of the polarization tensor within the su- 
perfluid phase were carried out at the one-loop approximation. This treatment was 
recently shown to be inadequate for the S-wave superfluid in neutron-star crusts. 

A many-body framework that is consistent with the sum rules for the polariza- 
tion tensor, in particular with the / sum rule 

i^jcLiu Imn^(q, w) = 0, (10) 

is provided by the random-phase resummation of the particle-hole diagrams in the 
superfluid matter. Because of the Nambu-Gorkov extension of the number of pos- 
sible propagators in the superfluid phase, which now include both normal (G) and 
anomalous (F) ones, at least three topologically different vertices are involved, which 
obey (schematically) the following equations 

f 1 = To + viGTiG + FT3G + GT2F + FT4F), (11) 
f 2 = v{GT2G^ + FT^G^ + GTiF + FTsF), (12) 
fs = viG'TsG + FTiG + G^TiF + #r2#), (13) 



lim 
q 
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kp [fm-i] 


m* /m 


A [MeV] 


Tc [MeV] 


i?(0.5) 


i?(0.8) 


i?(0.9) 


0.8 


0.97 


3.15 


1.78 


0.014 


1.0 


6.5 


1.6 


0.84 


0.57 


0.38 


0.022 


1.0 


7.4 



Note: Quoted are the wave number of neutrons, their effective mass, the gap and critical temper- 
ature, and the ratio i? as a function of reduced temperature T/Tc. 

V being the scalar interaction in the particle-hole channel. The fourth integral equa- 
tion for the vertex r4 follows upon interchanging particle and hole propagators in 
Eq. (Ilip . The full polarization tensor is the sum of the contributions shown in Fig. [3l 
It can be expressed through "rotated" polarization functions A, B, and C as 

A{q)C{q) + B{q)2 
Ciq)-v^[Aiq)Ciq)+Biq)2] 



with A = 2A2 loiq) - A2^, lAiq), B = -uA Io{q), and C((?) = -(^72) /o('?) + 
Ic{q)j where — g^/2m is the nucleon recoil. (The integrals Iq, Ic, and I a can 
be found in Rcf. 14.) It is now manifest that (q) = when q = 0. Thus, the 
leading order contribution to the polarization tensor appears at O(g^) and is linear 
in ^q. Since the neutrinos are thermal, with energies u; |q| T, the polarization 
tensor is suppressed by a factor T/m, which is of order 5 x 10^'^. 

The emissivity of the pair-breaking process can be compared to that of the 
modified bremsstrahlung (MB) process n + n^n + n + v + D, which is suppressed 
by roughly a factor exp(— 2A/T) in the superfluid phase. Thus, the ratio of the 
neutrino loss rate through MB to that from the pair-breaking process, as computed 
by Friman and Maxwell, is 

where gA and cy are the weak axial and vector coupling constants, to* and m,r 
are the neutron and pion masses, and F ~ 0.6 (defined in Ref. 11). The factor 
K = 0.2 accounts for the correction to the one-pion-exchange rate due to the full 
resummation of ladder series in neutron matter. The pair-breaking process domi- 
nates the MB process for temperatures below O.STc, where it is most efficient. This 
is illustrated in the table above. The comparison made here should be taken with 
caution, since the exponential suppression of the MB rate is not accurate (within a 
factor of a few) at temperatures close to the critical temperature. Nevertheless, one 
may safely conclude that the vector-current pair-breaking process is competitive 
with the modified pair-bremsstrahlung process in the relevant temperature domain 
T/Tc e [0.2-1]. 



5. Closing remarks 

Subnuclear matter at finite temperatures offers a fascinating arena for the develop- 
ment of many-body theory. Since the interactions invoved are well constrained by 
experiment, the entire complexity arises from the many-body correlations. As shown 
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in the examples chosen, subnuclear matter may exhibit a range of sahent many- 
body phenomena, such as Bose-Einstein condensation of alpha particles, BEC-BCS 
crossover in the deuteron channel, many-body extinction of bound states with in- 
creasing degeneracy, and non-trivial and quantitatively important modifications of 
the weak interaction rates due to many-body effects. 
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